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ABSTRACT 


M54 is a prototype for a globular cluster embedded in a dark matter halo. Gaia EDR3 photometry 
and proper motions separate the old, metal-poor stars from the more metal rich and younger dwarf 
galaxy stars. The metal poor stars dominate the inner 50 pc, with a velocity dispersion profile that 
declines to a minimum around 30 pc then rises back to nearly the central velocity dispersion, as expected 
for a globular cluster at the center of a CDM cosmology dark matter halo. The Jeans equation mass 
analysis of the three separate stellar populations gives consistent masses that rise approximately linearly 
with radius to 1 kpc. These data are compatible with an infalling CDM dark matter halo reduced to 
c 3 x 105 M at the 50 kpc apocenter 2.3 Gyr ago, with a central globular cluster surrounded by the 
remnant of a dwarf galaxy. Tides gradually remove material beyond 1 kpc but have little effect on the 
stars and dark matter within 300 pc of the center. M54 appears to be a "transitional" system between 
globular clusters with and without local dark halos whose evolution within the galaxy depends on the 


time of accretion and orbital pericenter. 


1. INTRODUCTION 


'The relationship between the very old globular clus- 
ters and sub-galactic mass cosmological dark matter ha- 
los (Doroshkevich et al. 1967; Peebles & Dicke 1968; 
Peebles 1984), is important to understand the earliest 
phases of star and galaxy formation. The observed in- 
ternal kinematics of globular cluster stars over a range 
of stellar masses combined with point mass n-body mod- 
eling of individual clusters finds that the mass den- 
sity of the clusters is consistent with the stellar mass 
of an old, evolved, stellar population inside 2-3 half- 
mass radii (Zonoozi et al. 2011; Baumgardt & Hilker 
2018; Ebrahimi et al. 2020; Vasiliev & Baumgardt 2021). 
At larger radii the declining density of cluster stars 
makes comtamination by foreground and background 
stars an increasing problem for velocity dispersion mea- 
surements. However, the growing number of radial ve- 
locity measurements and the increasing time baseline of 
the Gaia data (Gaia Collaboration et al. 2016) is allow- 
ing the kinematics of stars at larger radii to be examined. 

'There is evidence that some, but certainly not all, 
halo globular clusters have velocity dispersion profiles 
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at large cluster-centric radii above the predictions from 
the stellar mass of the cluster alone. Comparisons of the 
observed profiles to detailed n-body velocity dispersion 
profiles find that ~10% of currently measured clusters 
have velocity dispersion profiles above the values pre- 
dicted from the stars (Baumgardt & Hilker 2018; Bian- 
chini et al. 2019; Wan et al. 2021). For a set of 19 halo 
clusters measured to at least 5 half mass radii about 10% 
have a rising velocity dispersion profile and another 30% 
are approximately flat (Carlberg & Grillmair 2021). The 
current measurements are still somewhat noisy but will 
improve significantly with future Gaia releases. 

There are a number of possible explainations for some 
or all of the excess velocity dispersion, including single 
epoch observations of binary stars (Wan et al. 2021), the 
complexity of orbits between the cluster and the tidal 
surface (Fukushige & Heggie 2000; Claydon et al. 2017), 
particularly for highly elliptical orbits near the galactic 
center (Carlberg & Grillmair 2021), stripped galactic 
nuclei (Kuzma et al. 2018; Wirth et al. 2020), and dark 
matter halos (Mashchenko & Sills 2005a,b; Peñarrubia 
et al. 2017; Boldrini et al. 2020; Vitral & Boldrini 2021; 
Carlberg & Keating 2022). Dynamical measurements of 
potential progenitor systems can be used to test these 
ideas. 

The M54/NGC6715 cluster is at the center of the 
Sagittarius (Sgr) dwarf galaxy (Ibata et al. 1994; Bel- 
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lazzini et al. 2008). M54, with a mass of 1.4 — 1.8 x 
10°Mo, is one of the most massive globular clusters 
within the Milky Way halo. M54 satisfies the defini- 
tion of nuclear star cluster (NSC), one of four within 
smaller galaxies that have been identified as having ac- 
creted into the Milky Way (Pfeffer et al. 2021). As a nu- 
clear star cluster M54 is around the median NSC mass 
and one of the most metal poor known (Neumayer et al. 
2020). The stars in M54 are predominantly old and 
metal poor (Bellazzini et al. 2008; Alfaro-Cuello et al. 
2019) as opposed to several nearby NSCs whose central 
regions are predominantly young stars (Hannah et al. 
2021). M54 has no compelling evidence for a central 
black hole (Ibata et al. 2009; Wrobel et al. 2011). The 
Sagittarius dwarf luminosity, Ly = 1.8 x 109 Lo (Ma- 
teo 1998) is a factor of 10 less luminous than expected, 
relative to the mean relation for NSC hosts (Neumayer 
et al. 2020). The Sgr dwarf has prominent tidal tails 
(Ibata et al. 1997; Belokurov et al. 2006; Ibata et al. 
2020) which help constrain the recent orbital history of 
the system. Spectroscopic observations of the line of 
sight velocity profile of M54 find that the old, metal 
poor, blue stars of M54 have a velocity dispersion pro- 
file that declines from a central value of approximately 
10 kms-! to a minimum of about 5 kms-! at 22 pc 
and then rises back up to nearly 10 kms~! again at 55 
pc, where it becomes essentially the same as the velocity 
dispersion of the younger, more metal rich, redder stars 
of the dwarf galaxy surrounding the cluster (Bellazzini 
et al. 2008; Alfaro-Cuello et al. 2020). 

This paper uses Gaia data to separate the stars within 
about 4 kpc of M54 into blue, red and an intermedi- 
ate color (“green”) populations. The proper motions 
give the projected radial and tangential velocities of the 
stars. The resulting radial velocity dispersion profiles 
and the degree of radial anisotropy of each distinct stel- 
lar population are used with the Jeans equation to mea- 
sure the total mass profile with radius. The results are 
compared to n-body models to constrain the dark mass 
near and into the star cluster and within the surround- 
ing dwarf galaxy. The models are also used to assess 
the long term evolution of the M54 star cluster-dark 
halo systems and its dependence on the system's orbital 
pericenter. 


2. DATA 


We use the Early Data Release 3 (Gaia Collabora- 
tion et al. 2021) of the Gaia mission (Gaia Collabora- 
tion et al. 2016) to select Bp — Rp colors, G magni- 
tudes, parallaxes, proper motions, and their estimated 
errors, for stars within 10 degrees of the globular cluster 
M54/NGC6715, or 4.587 kpc at our adopted distance 


of 26.28 kpc (Baumgardt & Vasiliev 2021). The Gaia 
photometry is corrected for extinction using the red- 
dening maps of (Schlegel et al. 1998), themselves cor- 
rected using the prescription of Schlafly & Finkbeiner 
(2011), and using the Gaia DR2 coefficients derived by 
Gaia Collaboration et al. (2018). A Go, (GBP — Gre)o 
color-magnitude locus for the red giant branch (RGB), 
subgiant branch, and main sequence is constructed us- 
ing stars between one and three arcminutes from the 
cluster center. A theoretical isochrone with [Fe/H] — 
-1.59 (Harris 1996) generated in the Gaia passbands 
at http://stev.oapd.inaf.it/cgi-bin/cmd (Girardi et al. 
2004) is used as a guide at fainter magnitudes. This 
isochrone is not used directly due to significant and 
variable offsets between the observed and theoretical 
giant branches. Color offsets from this sampled color- 
magnitude locus are then computed using the tabulated 
photometric uncertainties in Gaia EDR3. The proper 
motion uncertainties lead to undertainties in each veloc- 
ity component of 5 and 10 kms~! at Go=15.5 and 17.0 
mag, respectively. The data are limited to an extinction 
corrected brightness of Go=17.5 magnitudes. Hereafter 
we drop the 0 subscript for the extinction corrected mag- 
nitudes and colors. 

There are several distinct age and metallicity stel- 
lar populations within the M54 cluster and Sagittarius 
dwarf (Bellazzini et al. 2008; Ibata et al. 2009; Vasiliev 
& Belokurov 2020). Detailed spectroscopic studies find 
that the cluster stars are predominantly an old metal 
poor population with a continuous spread to a younger 
metal rich population (Alfaro-Cuello et al. 2019). Along 
the lines of previous studies (Bellazzini et al. 2008; del 
Pino et al. 2021) we use the Gaia photometry to sep- 
arate the stars into three groups: the blue population, 
which is old and metal poor, a red population which is 
young and metal rich, and a population with intermedi- 
ate colors that we designate “green”. The blue stars 
are those having a color within 1 sigma of the clus- 
ter sequence. The red stars are defined as those be- 
tween Bp — Rp > 1.45 + 0.15(16 — G) and Bp — Rp < 
1.57 + 0.15(16 — G) for 18 > G > 16, and for G « 16 
colors between Bp — Rp > 1.45 + 0.5(16 — G) and 
Bp — Rp < 1.57 + 0.5(16 — G). The intermediate color, 
"green", stars are those that are at least 10 sigma to 
the red of the cluster sequence and bluer than the blue 
edge of the red sequence. The resulting color-magnitude 
diagram is shown in Figure 1. The sample shown in 
Figure 1 includes cuts on velocity relative to the cluster 
center of mass, as discussed below. 

'The M54 cluster is located close to the Galactic plane, 
at Galactic latitude of -14.1? and longitude 5.6? result- 
ing in a large and variable foreground of Galactic stars. 
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Figure 1. The extinction corrected color-magnitude dia- 
gram for the velocity selected star sub-samples. 
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Figure 2. The unit area normalized distribution of relative 
stellar velocities parallel to the galactic plane for the three 
color sub-samples. The top panel is for stars within 1 kpc of 
the cluster, the bottom for 4 kpc. The velocities perpendic- 
ular to the plane have similar distributions. 


The sample also has systematic errors with position (del 
Pino et al. 2021). The parallax errors are sufficiently 
large for most of the stars that a parallax limit is not 
helpful. Velocities are calculated assuming that the stars 
are at the distance of the cluster and placed in the center 
of mass frame of the cluster. Velocity cuts are important 
to reduce both the background and systematic errors. 

The radial velocities in the plane of the sky are projec- 
tion corrected, using the center of mass velocity of M54, 
143.06 kms! (Baumgardt et al. 2019) at a consensus 
distance of 26.28 kpc (Baumgardt & Vasiliev 2021). We 
will designate the two coordinates in the plane of the sky 
Y and Z, which are close but not quite identical to the 
Galactic Y and Z. Relative to the Galactic frame, our 
Z coordinate is tilted 5.5? backward from us and our 
Y coordinate is rotated 14? counterclockwise, as seen 
from above. The distribution of the velocities parallel 
to the Galactic plane, relative to the cluster, is shown 
in Figure 2 for stars within 1 and 4 kpc of the cluster 
(top and bottom panels, respectively). There is a clear 
excess of stars above 30 kms! at large radii. The ve- 
locities perpendicular to the plane show a similar excess 
of blue stars with velocities larger than about 30 km s7! 
beyond 1 kpc, most of the stars being those close to the 
plane. Gaia data has various sky location selection ef- 
fects in the region of M54 (see del Pino et al. (2021) for a 
discussion of the DR2 data). To produce an acceptably 
uniform sample in the M54 region we select stars that 
have a projected velocity relative to the cluster of 30 
km s^! or less. The distribution of the selected stars on 
the sky is shown in Figure 3. While Galactic plane stars 
intrude at the top of the figure, background contamina- 
tion has been greatly reduced. Using a larger velocity 
cut of 50 kms~! makes little difference for densities but 
boosts velocity dispersions at all radii, with a strong rise 
beyond 300 pc. We adopt the 30 kms! limit for our 
measurements. 


3. M54-SGR DENSITY AND VELOCITY 
DISPERSION PROFILES 


'The density and velocity dispersion profiles of the star 
cluster are measured as projected azimuthally averaged 
quantities, which is appropriate for the M54/NGC6715 
as it has an ellipticity of 0.06 (Harris 1996). Beyond 
100 pc the dwarf galaxy stars are in a significantly ellip- 
tical distribution, as shown in Figure 3. The Sagittar- 
ius dwarf galaxy is modeled as a prolate triaxial dwarf 
galaxy with its axis of rotation tilted approximately 40- 
60? from the sky plane (del Pino et al. 2021) with gra- 
dients of both density and velocity dispersion so shallow 
that radial averaging is acceptable. The radially av- 
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Figure 3. The sky distribution of the M54/Sgr |vxy| < 
30kms-! stars plotted in our longitude and latitude coor- 
dinates converted to kiloparsecs at the adopted distance of 
M54. The Galactic mid-plane is at y—6.46 kpc in these co- 
ordinates. The inner +0.1 kpc region, marked with a square 
at the center, is shown in the inset. At the distance of M54, 
1 kpc is 2.18? or 7.64 pc per arc-minute. Crowding causes 
an increasingly large under-sampling in the central region, 
where the M54 stars dominate those of the Sgr dwarf. 


eraged surface densities of the three color samples are 
shown in Figure 4. 

The red and green stars dominate the surface density 
beyond 50 pc in the Sgr dwarf, with blue stars com- 
prising about 5% of the combined surface density be- 
yond 100 pc. The central part of the field is so crowded 
that the density profile of M54 is severely undersampled 
inside 10 pc (Pancino et al. 2017), nearly twice M54’s 
half-mass radius of 5.6 pc. The green and red subsam- 
ples show only a weak rise towards the center in these 
data. Better resolved data for the young and metal rich 
red stars finds that their density rises within the cluster 
approximately in proportion to the blue stars, but with 
about 10% of the surface density (Bellazzini et al. 2008). 
Sgr is unusual in having a dominant old, metal poor pop- 
ulation, whereas young, metal rich stars often dominate 
the central regions of well resolved nearby nuclear star 
clusters (Hannah et al. 2021). 

The 2D velocity dispersion in the plane of the sky is 
oan = y o2 + o2, where o, and c; are the velocity dis- 
persion measured in the plane of the sky relative to the 
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Figure 4. The azimuthally averaged surface density versus 
projected radius for the color selected stars with absolute ve- 
locities less than 30 kms~!. Image crowding greatly reduces 
completeness inside approximately 20 pc. The central M54 
cluster has a 3D half mass radius of 5.6 pc. 


cluster center. The velocity dispersions have the mean 
proper motion velocity error in each bin, which range 
from about 8 to 12 kms~!, subtracted in quadrature. 
The resulting projected radial velocity dispersion, c; (r) 
and tangential velocity dispersion, o, of the three color 
subsamples are shown in Figure 5. The tangential veloc- 
ity dispersion in the range of 50 to 100 pc is consistent 
with zero, with less than 10 data points per bin in this 
radial range making the measurement of a small velocity 
dispersion difficult. The combined 2D velocity disper- 
sion is shown in Figure 8 compared with the n-body 
velocity dispersions. 

'The red and green stars, primarily in the dwarf galaxy, 
have similar velocity dispersion profiles, with a shallow 
decline toward the center inside 1000 parsecs. The blue 
stars, primarily in the central star cluster, have a dis- 
tinctive velocity profile, with a minimum at 30 parsecs 
and beyond that a rise back up to nearly the central 
value. The velocity dispersion dip in the blue stars and 
the flat profile of the red stars was noted in the line- 
of-sight radial velocity data of Bellazzini et al. (2008); 
the data here confirm that behavior using the orthogo- 
nal components provided by the Gaia proper motions. 
The lower velocity dispersion of the red stars relative to 
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Figure 5. The error corrected projected radial (top panel) 
and tangential (bottom panel) velocity dispersion profile of 
M54/Sgr for the blue, green, and red subsamples. The 2D 
velocity anisotropy, Bap = 1 — 07/02, evaluated in wider 
radial bins, is shown in the inset using the same line colors. 
The tangential velocity dispersion between 50 and 100 pc, 
where there are less than 10 data points per smaller bin, is 
consistent with zero 


the intermediate green stars was seen in the line-of-sight 
velocity analysis of Alfaro-Cuello et al. (2020). 

The inset in Figure 5 shows the 2D velocity anisotropy, 
Bop = 1 — o2/o2, in somewhat wider radial bins to in- 
crease the signal-to-noise ratio. The velocity ellipsoid of 
the blue stars becomes nearly completely radial around 
80 pc. At larger radii all three color subsamples have 
nearly isotropic velocity ellipsoids. The blue stars in the 
central 50 pc are kinematically distinct from stars of all 
colors at larger radii, which belong to the Sagittarius 
dwarf galaxy and have ages ranging from the last 0.5 
Gyr to a Hubble age (Alfaro-Cuello et al. 2019). 


4. TOTAL MASS PROFILE WITHIN 1 KPC 


The radial Jeans equation measures the total mass 
interior to radius r, 


ro? fdinv dino? 


M = 
Si G |dlnr S dlnr 


+26(r)}, (1) 


where v(r) is the density of a tracer population, o;.(r) is 
its radial velocity dispersion, and (r) is its 3D velocity 
dispersion anisotropy. The data displayed in Figures 4 
and 5 are used in the Jeans equation to estimate total 
interior masses shown in Figure 6. The Jeans equation 
assumes that the stars are in equilibrium and we are 
making the assumption of a spherical underlying po- 
tential. These assumptions are examined with models 
below. 

The power-law slope, s, = dlnv/dlnr, of the 3D stel- 
lar density profiles, v(r), is estimated from the slope of 
the surface density profile of the population, using the 
power law relation that the 3D slope is the 2D slope plus 
one. The 2D f values from the Figure 5 inset are close 
to the 3D values, 63p = 1 — (o2, + o2 )/202, for our 
assumed spherically symmetric system in which the two 
components, ou and ou, of the 3D tangential velocity 
dispersion, are equal. For the two special cases of the 
velocity ellipsoid being isotropic or completely radial, 
Bop = Bap, with values of 0 and 1, respectively. 

Figure 6 shows the results of the Jeans mass analysis. 
An important outcome is that the three separate star 
samples give consistent results at 30 pc and beyond. 
'The mass inside 30 pc is estimated to be 1.4, 1.9 and 
1.1x109 M for blue, green and red stars, respectively. 
'The 30 pc mass is consistent with the stellar dynamical 
mass of the star cluster alone, though it does not rule 
out a small fraction of dark matter. Between 30 and 95 
pc the total interior mass rises from the average value 
of 1.5 x 108 to 3.5 x 109 Mo, a factor 2.3. There is little 
increase in the stellar mass in this region, so the bulk of 
this mass increase must be made up of dark matter. The 
mass rises another factor of about 17 from 95 to 950 pc 
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Figure 6. Mass estimates from the Jeans equation for the 
three independent color subsamples. The mass points par- 
tially overlap at the largest radii. The star symbol indicates 
the Baumgardt & Hilker (2018) stellar mass of the cluster 
which is placed at 3 half-mass radii. 


to 89 x 109 Mo. Again, this must be made up mostly of 
dark matter, though the stars of the Sagittarius dwarf 
(with about 7096 of the luminosity Ly = 1.8 x 10" Lo 
within 950 pc) will contribute about 15% of the mass 
for a stellar M/Ly ~ 1. 


5. COMPARISON WITH N-BODY MODELS 


Our model for the M54-Sgr system begins when the 
system has settled into its current orbit, likely hav- 
ing shed the outer part of a more massive dark mat- 
ter halo (Lokas et al. 2010). A King model star clus- 
ter (King 1966) is inserted into a Hernquist model dark 
halo (Hernquist 1990). The two populations equilibrate 
in 0.02 Gyr. The combined star cluster-dark halo sys- 
tem is placed at the apocenter of M54’s orbit and inte- 
grated backward in time in a MW2014 potential (Bovy 
2015) with the bulge component modeled as a 0.5 kpc 
Plummer sphere. The large mass-to-light ratio of the 
Sagittarius dwarf, 22 Mo/Lo (Mateo 1998) or the 15% 
within 1 kpc we found in our Jeans mass analysis means 
that the stellar mass of the dwarf galaxy contributes 
little to the overall potential and can be ignored. The 
initial system is spherically symmetric. We select an in- 
ner, flattened, tilted, subset of dark matter particles on 
relatively circular orbits as proxies for the dwarf galaxy 
disk stars to follow their evolution with time. 


5.1. The N-body Model Star Cluster 
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Figure 7. The n-body model cluster stars (blue) and one 
of every 50 dark matter particles in the 3 x 10°Mo halo 
simulation. 


We adopt a star cluster mass of 1.78 x 109 Mc and 
a halfmass radius of 5.6 pc for M54 (Baumgardt 
& Hilker 2018) (https://people.smp.uq.edu.au/ Holger- 
Baumgardt/globular). The star cluster along with the 
surrounding dark matter is integrated with Gadget-4 
(Springel et al. 2021) modified to include a collisional 
star cluster (Carlberg & Keating 2021). Here we have 
included stellar “kicks” (Meiron et al. 2021) to mimic 
binary encounters, modeled as the addition of Gaussian 
random velocities of o \/3/2 added to each component of 
the velocity for a random fraction of 0.001 of the stars 
within half of the virial radius, where c is the 1D ve- 
locity dispersion. The cluster expands with time, so we 
start the simulation with half mass radius of 5.0 pc. The 
cluster expands to a half-mass radius of 6.6 pc at 2.29 
Gyr although this makes little difference to the dynam- 
ics beyond 10 pc. 

The dark halo is modeled as a Hernquist sphere (Hern- 
quist 1990) specified with a mass and scale radius with 
values motivated by the sub-galactic halos found in our 
cosmological simulations (Carlberg & Keating 2021). 
'The pre-infall mass of the entire Sagittarius galaxy to its 
virial radius may have been quite large, 1.6 x 1019 M. 
(Lokas et al. 2010) but we are interested in the mass 
remaining on the current orbit. The dark matter mass 
within 10 pc of the globular cluster center ( about two 
half-mass radii) is observationally constrained to be well 
below the stellar mass, otherwise the velocity dispersion 
of the cluster stars would be larger than observed. The 
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dark matter mass within 100 pc must rise to be at least 
the mass of the stellar cluster in order to have any dy- 
namical effect. We run models with Mp= 1 x 105 Mo 
and 3x 105 Mo, both a = 0.4 kpc, and Mp = 7x 108 Ms, 
a = 0.6 kpc. The star particles have mass 20 Mo and 
the dark matter particles 40 Mo with softenings of 2 and 
5 pc, respectively, which is sufficient to provide an accu- 
rate star cluster model at the half mass radius and be- 
yond. The current location and velocity of M54 (Baum- 
gardt et al. 2019) is used to integrate the orbit backward 
to the orbital apocenter around 2.29 Gyr ago (Vasiliev 
& Belokurov 2020). 

The star cluster is placed at the center of the dark 
halo. It is possible that the cluster formed at some dis- 
tance from the center of the dark halo and that dynami- 
cal friction subsequently dragged it inwards. Dynamical 
friction would have caused the massive, dense, M54 to 
spiral to the center of the dark matter halo with minimal 
tidal mass loss. For instance, a straightforward simula- 
tion shows that if started on a circular orbit at 0.5 kpc, 
the star cluster spirals into the center in less than 0.5 
Gyr, heating the central dark matter particles to create 
a 30 pc core in their distribution. 

'The star cluster and its dark halo adjust quickly to 
the new equilibrium with each other and their particle 
gravitational softening (2 and 5 pc, respectively). In the 
first 0.01 Gyr the star cluster half mass radius increases 
about 1096 and the dark halo density at 10 pc approx- 
imately doubles, with its velocity dispersion going up 
about 5096. 'The simulation responds dynamically on 
the radial orbit timescale, 0.86 Gyr. Successive peri- 
center passages heat the outer dark matter halo, which 
expands over the orbital period and is pulled away at 
the next pericenter passage. Tides and internal dynam- 
ical heating of stars moves about 1% of the cluster mass 
into the region between 30 and 200 pc. If there were no 
dark halo these stars would be unbound from the star 
cluster and would form a thin tidal stream of stars. In- 
ternal relaxation is very slow, as expected for M54 whose 
relaxation time has been estimated to be about 6 Gyr 
(Baumgardt & Hilker 2018). By comparison, our model 
has a relaxation time of 6.3 Gyr. 


5.2. Evolution of the Simulation 


'The model runs for 2.29 Gyr to reach a Galactic po- 
sition close to that of the observed cluster. The star 
particle positions and velocities are projected onto the 
sky and analyzed in the same was as the M54 cluster 
stars. The 2D velocity dispersion of M54's blue stars 
are compared to the star particles in the model in Fig- 
ure 8. The models predict a velocity minimum around 
30 pc and then a rise back to nearly the central veloc- 
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Figure 8. The projected 2D (top) and radial (bottom) ve- 
locity dispersion of M54’s stars compared to the same quan- 
tity measured for the star particles in the simulations that 
start with 1, 3, and 7x 10° dark halos, shown as black dashed, 
dotted and solid lines, respectively. 


ity dispersion at 100 pc. Our blue subsample shows the 
same trend. 

The model star cluster particles gradually extend to 
somewhat more than 100 pc. The M54/Sgr stars se- 
lected to be in the blue subsample extend over the en- 
tire surface of the dwarf galaxy, with cluster stars in the 
center and dwarf galaxy stars at larger radii. A kine- 
matic discriminant between cluster and Galaxy stars is 
the velocity anisotropy. Defining fp = 1 — o2/c2 the 
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Figure 9. The 2D velocity anisotropy, fap = 1— 02 / 02, for 
the color subsamples compared to the simulated star clus- 
ter's initial dark halo. The velocity dispersions have been 
corrected for the proper motion velocity errors. The values 
in the 1, 3, and 7x10? dark halos are shown as black dashed, 
dotted and solid lines, respectively. 


model values are compared to the color subsamples in 
Figure 9. The cluster stars in the model are expected to 
be radially anisotropic in the 30 to 100 pc range, as their 
orbits are generally perturbed only in the central region 
of the cluster and become highly elliptical. The green 
and red samples are essentially isotropic at all radii. 
The M54/Sgr data show that the blue sample is radially 
anisotropic in the range 30-100pc, becoming isotropic at 
larger radii. The change in velocity anisotropy with ra- 
dius is also consistent with the change in density profile 
shown in Figure 4, with stars inside of 100 pc largely be- 
longing to the cluster while stars beyond 100 pc largely 
belong to the dwarf galaxy. 


6. THE MEAN VELOCITY FIELD OF THE 
SAGITTARIUS DWARF 


'The Sagittarius Dwarf galaxy is being shredded in the 
Milky Way's gravitational field (Ibata et al. 1994) to 
produce a tidal stream wrapped around the galaxy at 
least twice (Ibata et al. 2020). Our interest is to use the 
velocity profile of the inner few kiloparsecs of Sagittarius 
stars to estimate the density profile of the dark matter 
that remains bound to Sagittarius and place some limits 
on the pre-accretion mass of the system and how it is 
related to the M54 cluster. The stars in this region 
have motions that transition from bound to flowing out 


into the stream. The velocity field of the Sagittarius 
Dwarf has been extensively studied to understand the 
shape of the stellar system, the total mass profile and 
the velocity at which stars are fed into the stellar stream 
(Law & Majewski 2010; Penarrubia et al. 2010; Lokas 
et al. 2010; Frinchaboy et al. 2012; Vasiliev & Belokurov 
2020; del Pino et al. 2021; Ramos et al. 2021). A point of 
controversy is the amount of rotation in the inner 0.5-2 
kpc of the system. 


6.1. The Sagittarius Galaxy Velocities 


Our color subsamples are separately binned into 0.4 
kpc square pixels in a 8 kpc frame centered on the M54 
cluster. The mean velocity fields on the plane of the sky 
(YZ Galactic coordinates relative to the cluster center) 
are shown in the upper panels of Figures 10, 11 and 12. 
The motions are modeled as the sum of a rotation and 
a shear, 


T 


Ut (r) — Urot JEFE (2) 
vr(r)=Arr (3) 
U(r, 0) = Asr|cos (0 — 0o), sin (0 — Beil, (4) 


which are projected into the x and y components of the 
velocity. The parameters of the model are A, the lin- 
ear expansion rate, Vrot the peak of the maximum of 
the rotation velocity, a the scale radius of the rotation 
velocity, A, the velocity shear, and ĝo the direction of 
the shear. The model is fit to the data using the Monte 
Carlo Markov Chain algorithm emcee (Foreman-Mackey 
et al. 2013). The logarithm of the likelihood function is 
—$ > [A (v)? + ^(v)?]/o2 + log (02), where A(v,) is 
the difference between the observed v, and the value 
that the model predicts. The c; is taken from the mea- 
sured 2D velocity dispersion which is essentially uniform 
beyond 0.2 kpc from the cluster, nearly independent of 
color. We set c to be constant at 20 kms-!. The pri- 
ors on the model quantities are set to be flat with such a 
large range that the results have no dependence on the 
range of the prior. 

'The resulting corner plots showing the fits are dis- 
played in the lower panels of Figures 10, 11 and 12. The 
model accounts for 44 to 6196 of the variance in the data 
overall so the fits are reasonable, with y?/dof values 
ranging from 1.1 to 2.1 assuming a random velocity of 20 
kms~! in each component. Nevertheless, the model does 
indicate that all three subsamples have similar mean mo- 
tions, with an expansion of A, = 1.6—2.0kms ! kpc t 
and a shear of A, = 2—2.5kms ! kpc !, at an angle of 
40-45? with respect to the Galactic plane. 


6.2. An Artificial Star Dwarf Galaxy 
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Figure 10. The green sample mean velocity field (top panel) 
and the corner plot of the fit to a rotation plus shear model. 
The model removes 73% of the variance and has a x?/dof = 
2.06 


To approximate the kinematics of the stars of the 
Sagittarius dwarf we select a thick disk-like set of parti- 
cles from the simulation at early times. We select from 
the simulation near the outset, at an age of 0.05 Gyr, 
identifying a subset of the dark matter particles around 
the star cluster which are in a flattened distribution, 
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Figure 11. The blue sample mean velocity field (top panel) 
and the corner plot of the fit to a rotation plus shear model. 
The model removes 55% of the variance. x?/dof = 1.10 


tilted to the line of sight. Sagittarius is a low surface 
brightness galaxy in which the stars provide relatively 
little mass compared to the dark matter (Ibata et al. 
1997; Mateo 1998). The model dwarf galaxy is con- 
structed 2.24 Gyr ago shortly after the system passed 
apocenter. The initial disk orientation is taken to be the 
same as the current orientation, with a unit axis [-0.13, 
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Figure 12. The red sample mean velocity field (top panel) 
and the corner plot of the fit to a rotation plus shear model. 
The model removes 66% of the variance. x?/dof = 2.02 


0.96, 0.26] (del Pino et al. 2021). The artificial stars are 
selected as dark matter particles within 4 kpc of the cen- 
ter in the plane of the disk, with a height above the disk 
less than 0.3 of their radius and a vertical velocity of 
15 kum e" or less. The 0.05 Gyr time mass profile of the 
1-3 and 7 x 105 M, halos are reasonably approximated 
with Hernquist spheres with scale radii of 0.35 and 0.55 
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Figure 13. The M54-Sgr stellar surface densities with ra- 
dius (blue, green, red), the n-body star cluster surface den- 
sity (blue dashed), the artificial dwarf galaxy star profiles 
(purple) and the dark matter surface density in units of 
10? Mepe~? (gray), where the 1, 3 and 7 x 105 Ms simu- 
lations have dashed, dotted and solid lines, respectively. 
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Figure 14. The mean 3D density inside radius r inferred 
from the Jeans equation in Figure 6 and the 3D dark mat- 
ter densities of the 3 and 7 x 10°Mo simulations as dotted 
and solid lines, respectively. The solid grey line is the lo- 
cal 3D density in the 7 x 105 M simulation. All simulation 
measurements are at 2.29 Gyr. 
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Figure 15. The starting and final radii of the artificial dwarf 
galaxy stars for the 3 x 10°Mo simulation. No particles 
are selected as stars inside 0.5 kpc, but the particles are on 
sufficiently elliptical orbits that over time the central hole 
fills. 
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Figure 16. The stream that the artificial dwarf galaxy star 
particles produce in the 3 x 105 Ma simulation. The Galactic 
center is located at (0,0) and the star cluster is located close 
to the observed position of M54, at Galactic coordinates [- 
14.7, 2.2,-6.2] kpc. The simulation began at the orbital apoc- 
enter 2.29 Gyr in the past. A dwarf galaxy rotating in the 
opposite direction does not produce a bifurcated stream. 


kpc, respectively, slightly smaller scale radii than the 
initial values. The dark matter halo mass model defines 
a circular velocity with radius and the energy and angu- 
lar momentum of circular orbits with radius. At 1 kpc 


10 40 60 80 


the circular velocity in the halo is 23 kms~!. Particles 
with energies between 0.8 and 1.2 of the circular veloc- 
ity energy and 0.8 to 1.2 of the circular velocity energy 
within Tmax of 4.0 kpc are selected at an initial time of 
0.05 Gyr as an artificial first approximation to the dis- 
tribution of dwarf galaxy stars. To create a flat density 
profile we select particles with rmin > 0.5 (1.0) kpc, for 
the 3 and 7 x 105 M; models, respectively. At the time of 
selection the particles are in a tilted, rotating disk of 4.0 
kpc extent and an aspect ratio of about 1:3. The mean 
rotation of the stars is counterclockwise around the an- 
gular momentum vector, but clockwise when viewed in 
the plane of the sky. 

'The projected surface density profile of the n-body 
stellar cluster (dashed blue line) is shown in Figure 13. 
The comparison with the density profile of the old, 
metal-poor stars of our Gaia subsample illustrates the 
tremendous undersampling of stars in the central region 
of the star cluster and how the stars of the cluster over- 
lap with the old, metal poor, stars of the dwarf galaxy 
beyond 50 pc. Three halo models, with initial masses 
of 1, 3 and 7x10* M5, were evolved along with the star 
cluster, with the projected surface density of the dark 
matter halos at the time that matches the current or- 
bital location are shown as the dashed, dotted and solid 
gray lines, respectively. Randomly selecting particles 
from the entire dark matter halo will not mimic the Sgr 
dwarf density profile, which requires that the fraction of 
dark matter particles identified as stars decline inward 
about a factor of 100 from 1 kpc to 10 pc. The simple 
identification procedure here leads to the surface densi- 
ties shown as purple lines in Figure 13. The procedure 
with an initial central hole in the selected stars leads to 
a somewhat too flat density profile which could be fixed 
by allowing some of the stars in the central region. The 
main outcome is that the dwarf galaxy stars are a de- 
creasing fraction of the dark matter towards its center. 

In Figure 14 the mean interior density, p(r) — 
M(r)/(4n/3r?) of the dark matter (without the star 
cluster mass) derived from the three independent color 
subsamples is compared to the three n-body dark mat- 
ter mean interior halo densities at the current location of 
the M54-Sgr system. Within this limited range of these 
three models the initial dark matter halo of 3 x 105 Mọ 
and scale radius 0.4 kpc gives results closest to the ob- 
served values. 

The distributions of the artificial dwarf galaxy star 
particles are shown in Figure 16 for the simulation with 
a total halo mass of 3 x 105 Mo. The Galactic XZ plane 
is shown since the orbit is nearly in this plane. The 
stream in the Mp = 7 x 105 M model is slightly wider, 
longer, and fuzzier. The lower mass model produces 
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Figure 17. The radial (top) and tangential (bottom) veloci- 
ties with radius for the M54/Sgr data and the artificial stars, 
selected with positive plane-of-the-sky tangential velocities. 
The central model cluster star particles are shown as dotted 
blue lines. The other lines are as in Figure 13. 


a significantly narrower stream. A similar plot, with 
the direction of the X axis reversed, is shown in Figure 
10 of Vasiliev & Belokurov (2020). Their simulations 
are based on dark halos with circular velocity curves 
that peak in the range 2-4 kpc, whereas here it peaks at 
around 0.4 kpc. 

It is interesting to note that the selection of positive 
tangential velocities (projected into the plane of the sky) 
leads to stream bifurcations, as seen in the Sagittarius 
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Figure 18. The mean velocity vectors in the plane of the 
sky of the artificial dwarf galaxy particles in the 3 x 108 Mo 
simulation. The plot is is comparable to velocity map panels 
of Figures 10, 11 and 12. The arrow in the upper left has a 
length of 30 kms" !. 


stream (Belokurov et al. 2006; Koposov et al. 2012). The 
same selection with negative velocities does not have 
prominent bifurcations, confirming the importance of ro- 
tation in the dwarf galaxy for their origin (Ramos et al. 
2021). 

One of the uses of the artificial dwarf galaxy stars is to 
identify effective tidal radii, within which stars remain 
bound to the dwarf and beyond which they eventually 
join the massive tidal stream. The radii of stars relative 
to the star cluster at 0.05 Gyr and 2.29 Gyr are displayed 
in Figure 15. The artificial star particles that join the 
stream have a clearly defined minimum radius, about 
0.5 kpc, for 3 x 109 M. halo. Beyond the minimum 
radius an increasing fraction of particles join the stream 
with radius, with almost all joining at 2 kpc for the 
3 x 109 M5 halo. The absence of a crisp tidal radius is a 
consequence of the elliptical Galactic orbit of the system 
and the eccentricity of the particle orbits. 

The velocities and positions of the artificial dwarf 
galaxy stars are projected onto the sky as seen from the 
location of the Sun in the Milky Way at time 2.29 Gyr for 
the two halo simulations. The velocities are measured in 
the same way as for the observational data. The result- 
ing, radially averaged velocity fields are shown in Fig- 
ure 17. Figure 17 shows the mean radial velocities in the 
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top panel and the corresponding model measurements. 
'The positive bump in the measured mean radial velocity 
of the blue stars has a large uncertainty but is also seen 
at a somewhat larger radius in the n-body model star 
cluster particles, where it is likely a result of the Galactic 
tidal field for the system which passed pericenter 0.04 
Gyr before the time being analyzed. Figure 17 shows 
the mean tangential velocities in the bottom panel. The 
3 x 105 M model has tangential velocities closer to the 
observed values over the 100-1000 pc range, leading us 
to favor it as the better velocity model inside 1 kpc. 

Beyond 1 kpc the strong shear seen in Sgr shows up as 
the outer tangential velocity decline in Figure 17, which 
is not present in the projection of any of the models. The 
artificial disk star velocities plotted over the face of sys- 
tem shown in Figure 18 are similar in the inner 1-2 kpc 
but not a good fit to the velocities of Figures 10, 11 and 
12 beyond 2 kpc. The model system has a smaller shear 
velocity from top to bottom than observed in Sagittar- 
ius. 

'The model system has a large velocity gradient along 
the line-of-sight velocity at the current orbital phase as 
can be seen in Figure 16, where the segment of the 
stream below M54 is moving nearly horizontally away 
from the galactic center, whereas the part of the system 
at the top, closer to the plane, is turning upward with 
a reduced line of sight velocity. To increase the shear 
in the model one solution would be to rotate the orbit 
about 5?, to project more of the strong velocity shear 
into the plane of the sky. Another possibility is that the 
cylindrically symmetric potential used here does not al- 
low for the Galactic bar, although the asymmetry of the 
potential due to the bar will be small at a pericenter 
distance of 14 kpc. 


7. ALTERNATIVE FUTURES EVOLUTION 


Our radial velocity dispersion measurements, their 
population density slopes, and @ values lead to the mass 
estimates shown in Figure 6. Our Jeans equation mass 
analysis finds 1.4, 1.9 and 1.1x 109 M within 30 pc, for 
the blue, green and red samples, respectively. The Jeans 
masses are essentially consistent with the cluster's stel- 
lar mass of 1.41 — 1.78 x 109 Mc (Baumgardt & Hilker 
2018). At 95 pc the Jeans analysis masses are 6.3, 5.9 
and 3.7 x 109 Mo, and at 300 pc the masses are 15.9, 18.2, 
and 13.8x109 M, for the blue, green and red samples, 
respectively. The implied rise of mass with radius is r!:2° 
between 30 and 95 pc and r!-°8 between 95 and 300 pc, 
with a total increase in mass a factor of 10.9 between 
30 and 300 pc. The mean density of the dark matter is 
close to what is expected for a CDM dark matter halo 
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Figure 19. Top: The 3D velocity dispersion at 10 Gyr for 
M54's orbit (solid) and the orbit started at apocenter with 
1/2 (dashed) and 1/4 (dotted) of the angular momentum. 
'The black line is for the dark matter, the blue for the cluster 
stars. Bottom: the density profiles for the same models with 
the same line designation. 


of around 3 x 105 M; initiated at the apocenter of the 
frictionless M54 orbit 2.29 Gyr ago. 

'The mass and extent of the remnant globular cluster 
dark matter halo depends on how long ago the system 
fell into the Milky Way and how close the orbit comes 
to the center of the galaxy, where the tidal field dynam- 
ically heats the system's particles and pulls them away 
in a tidal stream. Over the 2.29 Gyr of evolution the 
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models lose about 1696 of the mass within 1 kpc to tidal 
fields. It is not unlikely that systems like Sagittarius 
with a massive globular cluster near their centers fell 
into the Milky Way much earlier when the bulk of the 
dark matter mass was built up, which will lead to more 
tidal mass loss. Simply running the simulation to 10 
Gyr we find that the 3 x 105 M. halo loses 38% of its 
mass inside 1 kpc, but only 996 inside 0.1 kpc. Had the 
orbit come closer to the Galactic center even more dark 
matter mass would have been lost. Two further sim- 
ulations start at the same apocenter as used for M54, 
but reduce the initial orbital velocities so that the sys- 
tem has 1/2 and 1/4 of the angular momentum. The 
simulations use the same static Galactic potential with 
no allowance for dynamical friction. The resulting star 
and dark matter velocity and density dispersion profiles 
are shown in Figure 19. As expected, smaller pericenter 
passages remove more dark matter more quickly. The 
velocity dispersion profile of the system on the 1/4 an- 
gular momentum orbit is always declining to 100 pc. 
Tidal fields remove 9096 of the mass inside 1 kpc within 
the first 2 Gyr. If gas were present, it is likely that ram 
pressure would remove much of the gas as well (Tepper- 
García & Bland-Hawthorn 2018). Therefore, it appears 
that if M54 had fallen in earlier on a smaller pericen- 
ter orbit it would have formed few if any surrounding 
dwarf galaxy stars and lost so much of its dark halo that 
it would be observed today as a completely old globu- 
lar cluster with no detectable dark halo. On the other 
hand, if M54's and the surrounding halo was still in the 
galactic outskirts, it could have grown into a substantial 
dwarf galaxy and the central star cluster could well have 
grown into a much larger nuclear star cluster, possibly 
even into one containing a detectable black hole. In that 
sense, M54 can be considered a "transition" object. 


8. DISCUSSION AND CONCLUSIONS 


Our photometric selection using isochrones adjusted 
to the Gaia data separates the stars into blue, red (metal 
rich and younger) and an intermediate, “green” popula- 
tion. Each population is kinematically and spatially dis- 
tinct, although the red and green populations are sim- 
ilar. A Jeans equation mass analysis using the three 
populations finds that the mass inside 95 pc is approx- 
imately 3.5 times the stellar mass of the M54 star clus- 
ter. Beyond 100 pc the mass continues to grow approx- 


imately linearly with distance. The overall rise of mass 
from 30 to 1000 pc is nearly linear, M œ zl). implying 
an approximately r~? density profile, with a core radius 
no larger than 30-50 pc. N-body simulations that start 
2.29 Gyr in the past with a star cluster in a Hernquist 
sphere of masses 1 — 3 x 105 M; and a scale radius of 
400 pc provide a good match to most of the properties 
of the system. The system may have been much more 
massive when it first encountered the Milky Way, before 
it spiraled into the current orbit. 

Most of the Milky Way halo globular clusters have no 
significant dark matter detected in their outskirts. How- 
ever, ~10% of the accessible halo clusters have been 
tentatively identified as being consistent with having 
a local dark matter halo (Carlberg & Grillmair 2021). 
Simulations show that approximately 25-3596 of clusters 
started in the centers of local halos retain significant 
dark matter as the Galactic halo assembles (Carlberg & 
Keating 2021). The M54-Sagittarius system appears to 
be an intermediary system that would have become a 
standard globular cluster had it fallen in earlier and its 
orbit taken it into regions of somewhat stronger tides. 
On the other hand, if the system had stayed in the out- 
skirts the surrounding dwarf galaxy would likely have 
continued to form stars to become a more prominent 
dwarf galaxy today. 

The Sagittarius galaxy and its stream contains 8 other 
high confidence globular clusters (Da Costa & Arman- 
droff 1995; Bellazzini et al. 2020) including the distant, 
high mass NGC 2419 cluster (Ibata et al. 2013). Four 
of the eight are old, metal poor clusters, the other 4 
are more metal rich and younger, along the lines of the 
disk and halo clusters of the Milky Way (Searle & Zinn 
1978). An additional 10-20 lower luminosity clusters are 
proposed members of Sagittarius (Minniti et al. 2021). 
If the Milky Way is assembled out of systems similar to 
Sagittarius, then its one cluster (M54) in a dark halo of 
the nine luminous clusters somewhat weakly bolsters the 
idea that roughly 1096 of disk galaxy clusters retain pre- 
galactic dark matter halos. If the M54/Sagittarius sys- 
tem is a late-time accretion event that is representative 
of earlier ones (Massari et al. 2019; Malhan et al. 2022), 
then the globular clusters in those accretion events may 
still have one or two clusters with local dark matter ha- 
los orbiting freely. 


'This research was supported by NSERC of Canada. 
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